👉 Os testes de hipóteses em um contexto multivariado são muito mais complexos do que num contexto univariado.
Considerando a distribuição normal, temos
no caso univariado, uma média (\(\mu\)) e uma variância (\(\sigma^2\)), totalizando dois parâmetros,
no caso \(p\)-variado, \(p\) médias, \(p\) variâncias e \(\left( \begin{array}{c} p \\ 2 \end{array} \right)\) covariâncias, totalizando
\[p + p +\left( \begin{array}{c} p \\ 2 \end{array} \right) = \dfrac{1}{2} p(p+3)\]
parâmetros.
Por exemplo, para \(p = 10\), temos
\[\dfrac{1}{2} p(p+3) = \dfrac{10(10+3)}{2} = \dfrac{130}{2} = 65\]
parâmetros, para cada um dos quais uma hipótese poderia ser formulada.
Além disso, podemos estar interessados em testar hipóteses sobre subconjuntos desses parâmetros, ou ainda, sobre funções deles.
🤔 Por que utilizar testes multivariados ao invés de testes univariados?
\[\begin{eqnarray*} P(\text{pelo menos uma rejeição}) &=& 1 - P(\text{todos os 10 testes não rejeitarem } H_0) \\ &=& 1 - (0,95)^{10} = 0,40 \end{eqnarray*}\]
👉 Se as variáveis não forem independentes, teríamos \(0,05 < \alpha < 0,40\)
🤔 Por que utilizar testes multivariados ao invés de testes univariados?
Num primeiro momento, vamos revisar o caso univariado de se testar a média, \(\mu\), de uma população com \(\sigma^2\) conhecido.
Sejam as seguintes hipóteses:
\[H_0: \mu = \mu_0 \text{ vs. } H_a: \mu \neq \mu_0\]
Considere uma amostra aleatória \(X\) de tamanho \(n\), de forma que \(X \sim N(\mu, \sigma^2)\), com \(\sigma^2\) conhecido. Uma estatística de teste apropriada para este tipo de situação é
\[z = \dfrac{\bar{x}- \mu_0}{\dfrac{\sigma}{\sqrt{n}}} \stackrel{H_0}{\sim} N(0,1)\]
Em que \(\bar{x} = \dfrac{\displaystyle{\sum_{i=1}^n x_i}}{n}\). Para um nível \(\alpha\) de significância, rejeitamos \(H_0\) se \(|z| \geq z_{\frac{\alpha}{2}}\).
Equivalentemente, podemos usar \(z^2\) que segue uma distribuição \(\chi^2_1\) e rejeitamos a hipótese nula \(H_0\) se \(z^2 \geq \chi^2_1(\alpha)\).
Se \(n\) for grande, o Teorema Central do Limite garante que \(z\) é aproximadamente normal, mesmo que as observações não sejam normais.
Considere agora uma amostra aleatória \(\mathbf{x} = \{X_1, X_2, \cdots, X_n\}\) de uma população \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\), com vetor de médias \(\mathbf{\mu} = \{\mu_1, \mu_2, \cdots, \mu_p\}^t\) e matriz de variâncias e covariâncias
\[\mathbf{\Sigma} = \begin{bmatrix} \sigma_{11} & \sigma_{12} & \cdots & \sigma_{1p}\\ \sigma_{21} & \sigma_{22} & \cdots & \sigma_{2p}\\ \vdots & \vdots & \ddots & \vdots\\ \sigma_{p1} & \sigma_{p2} & \cdots & \sigma_{pp} \end{bmatrix}\]
Queremos testar as seguintes hipóteses:
\[H_0: \mathbf{\mu} = \mathbf{\mu}_0 \text{ vs. } H_a: \mathbf{\mu} \neq \mathbf{\mu}_0\]
Podemos generalizar a estatistica \(z^2\) do teste univariado da seguinte forma:
\[Z^2 = n(\mathbf{\bar{x}}-\mathbf{\mu}_0)^t\mathbf{\Sigma}^{-1}(\mathbf{\bar{x}}-\mathbf{\mu}_0)\]
Se \(H_0\) é verdadeira, \(Z^2 \stackrel{H_0}{\sim} \chi^2_p\) e, portanto, rejeitamos a hipótese nula \(H_0\) se \(Z^2 > \chi^2_p(\alpha)\)
Exemplo 01: Considere uma amostra aleatória de tamanho \(n = 20\), com vetor de médias \(\mathbf{\bar{x}} = (71,45;164,7)^t\), de uma população normal bivariada cuja matriz de variâncias e covariâncias seja dada por
\[\mathbf{\Sigma} = \begin{bmatrix} 20 & 100 \\ 100 & 1000\\ \end{bmatrix}\]
Suponha que desejássemos testar as seguintes hipóteses
\[H_0: \mathbf{\mu} = [70, 170]^t \text{ vs. } H_a: \mathbf{\mu} \neq [70, 170]^t\]
Temos que,
\[\begin{eqnarray*} Z^2 &=& n(\mathbf{\bar{x}}-\mathbf{\mu}_0)^t\mathbf{\Sigma}^{-1}(\mathbf{\bar{x}}-\mathbf{\mu}_0) \\ &=& 20 \times \left[ \begin{array}{c} 71,45 - 70 & 164,7 - 170 \end{array} \right]\begin{bmatrix} 20 & 100 \\ 100 & 1000\\ \end{bmatrix}^{-1}\left[ \begin{array}{c} 71,45 - 70 \\ 164,7 - 170 \end{array} \right] \\ &=& 20 \times \left[ \begin{array}{c} 1,45 & -5,30 \end{array} \right]\begin{bmatrix} 0,1 & -0,01 \\ -0,01 & 0,002\\ \end{bmatrix}\left[ \begin{array}{c} 1,45 \\ -5,30 \end{array} \right] = 8,4026 \end{eqnarray*} \]
Usando um nível de significância \(\alpha = 0,05\), temos que \(\chi^2_{(2;0,05)} = 5,991\) e, portanto, rejeitamos a hipótese \(H_0: \mathbf{\mu} = [70, 170]^t\), uma vez que \(Z^2 = 8,4026 > 5,991\).
n = 20
p = 2
alpha = 0.05
xbarra = c(71.45,164.7)
mi0 = c(70,170)
Sigma = matrix(c(20,100,100,1000), nrow = p, byrow = TRUE)
# H0: mu = [70, 170]^t - Estatística Z2
Z2 = n*t(xbarra-mi0)%*%solve(Sigma)%*%(xbarra-mi0)
# valor crítico da distribuição qui-quadrado
qui = qchisq(1-alpha, p)
# p-valor
pvalor = 1-pchisq(Z2, p)
Z2 > qui # Rejeita H0Considere o caso univariado de uma amostra aleatória de tamanho \(n\), \(X_1, X_2, \cdots, X_n\), da distribuição \(N(\mu,\sigma^2)\), em que ambos os parâmetros \(\mu\) e \(\sigma\) são desconhecidos.
Considere também, que desejamos testar se um determinado valor \(\mu_0\) é um possível valor para a média populacional \(\mu\). Do ponto de vista de teste de hipóteses, isso é equivalente a testar:
\[H_0: \mu = \mu_0 \hspace{0.5cm} \text{ vs } \hspace{0.5cm} H_a: \mu \neq \mu_0\]
Neste caso, uma estatística de teste apropriada é
\[t = \displaystyle{\frac{\bar{X} - \mu_0}{s/\sqrt{n}}}\]
em que \(\bar{X} = \dfrac{\displaystyle{\sum_{j=1}^{n}} X_{j}}{n}\) e \(s^{2} = \dfrac{\displaystyle{\sum_{j=1}^{n}}(X_{j}-\bar{X})^2}{n-1}\).
Esta estatística de teste tem uma distribuição t-Student com \(n − 1\) graus de liberdade. Rejeita-se \(H_0\), se o valor observado de \(|t|\) exceder um valor específico da distribuição t-Student com \(n − 1\) graus de liberdade.
Rejeitar \(H_{0}\) quando \(|t|\) é grande é equivalente a rejeitar \(H_{0}\) quando o valor de \[\begin{equation} t^{2}=\frac{(\bar{X}-\mu_{0})^{2}}{s^2/{n}}= n(\bar{X}-\mu_{0})(s^{2})^{-1}(\bar{X}-\mu_{0}) \label{eq1} \end{equation}\] for grande.
Uma vez que \(\bar{X}\) e \(s^{2}\) são observados, rejeitamos \(H_{0}\) em favor a \(H_{a}\) a um nível \(\alpha\) de significância, se: \[n(\bar{X}-\mu_{0})(s^{2})^{-1}(\bar{X}-\mu_{0})>t^{2}_{n-1,{\alpha}/{2}}\]
Ou, de forma equivalente,
\[n(\bar{X}-\mu_{0})(s^{2})^{-1}(\bar{X}-\mu_{0})>F_{1, n-1}(\alpha)\]
Se \(H_{0}\) não for rejeitada, concluímos que \(\mu_{0}\) é um possível valor para a média populacional \(\mu\).
Considere agora o problema de se determinar se um dado vetor \(\mathbf{\mu}_{0}\), \(p \times 1\), é um possível valor para o vetor de médias \(\mathbf{\mu}\) de uma distribuição normal \(p\)-variada.
Uma generalização natural da distância quadrática encontrada para o caso multivariado é dada por:
\[\begin{equation}T^{2}=(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\left ( \frac{\mathbf{S}}{n} \right )^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})=n(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})\end{equation}\]
em que
\[\bar{\mathbf{x}} = \frac{1}{n} \sum_{j=1}^{n}{\mathbf{X}_{j}}\]
\[\mathbf{S} = \frac{1}{n-1} \sum_{j=1}^{n}(\mathbf{X}_{j} -\bar{\mathbf{x}})(\mathbf{X}_j-\bar{\mathbf{x}})^{t}\]
\[\mathbf{\mu}_{0} = \left[ \begin{array}{c} \mu_{10}\\ \vdots \\ \mu_{p0}\\ \end{array} \right] \]
A estatística \(T^{2}\) é chamada de \(T^{2}\) de Hotelling.
Se o valor observado da estatística \(T^{2}\) é muito grande, isto é, se \(\bar{\mathbf{x}}\) está “muito longe” de \({\mathbf{\mu}}_{0}\), a hipótese \(H_{0}:{\mathbf{\mu}}={\mathbf{\mu}}_{0}\) é rejeitada.
Proposição 01: Sejam \(\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\) uma amostra aleatória de uma distribuição \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\) e \(T^2 = n(\bar{\mathbf{x}}-\mathbf{\mu}_{0})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-\mathbf{\mu}_{0})\). Então,
\[\displaystyle{\frac{n - p}{(n - 1)p}}T^2 \stackrel{H_0}{\sim} F_{p,n-p}\]
Dessa forma, um teste para as hipóteses
\[H_{0}:\mathbf{\mu} = \mathbf{\mu}_{0} \text{ vs. } H_{a} :\mathbf{\mu} \neq \mathbf{\mu}_{0}\]
pode ser formulado.
Ao nível \(\alpha\) de significância, rejeita-se \(H_{0}\) em favor de \(H_{a}\) se
\[T^{2}=n(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}}) > \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha)\]
de acordo com a proposição 01. \(F_{p,n-p} (\alpha)\) denota o \(100(1 − \alpha)\)-ésimo percentil superior de uma distribuição \(F[p; (n−p)]\).
Exemplo 02: Suponha que a matriz de dados para uma amostra aleatória de tamanho \(n = 3\) de uma população normal bivariada seja dada por
\[\mathbf{X} = \begin{bmatrix} 6 & 9 \\ 10 & 6\\ 8 & 3 \end{bmatrix}\]
Avalie o valor observado de \(T^2\) para \(\mathbf{\mu}_0 = [9,5]^t\). Qual a distribuição amostral de \(T^2\) neste caso?
Temos que
\[\mathbf{\bar{x}} = \begin{bmatrix} \bar{x}_1 \\ \bar{x}_2 \end{bmatrix} = \begin{bmatrix} \dfrac{6 + 10 + 8}{3} \\ \dfrac{9 + 6 + 3}{3} \end{bmatrix} = \begin{bmatrix} 8 \\ 6 \end{bmatrix}\] Além disso,
\[s_{11} = \dfrac{(6-8)^2 + (10-8)^2 + (8-8)^2}{2} = 4\]
\[s_{12} = \dfrac{(6-8)(9-6) + (10-8)(6-6) + (8-8)(3-6)}{2} = -3\]
\[s_{22} = \dfrac{(9-6)^2 + (6-6)^2 + (3-6)^2}{2} = 9\]
De forma que,
\[\mathbf{S} = \begin{bmatrix} 4 & -3 \\ -3 & 9 \end{bmatrix}\]
e,
\[\mathbf{S}^{-1} = \dfrac{1}{(4 \times 9) - ((-3) \times (-3))}\begin{bmatrix} 9 & 3 \\ 3 & 4 \end{bmatrix} = \begin{bmatrix} \frac{1}{3} & \frac{1}{9} \\ \frac{1}{9} & \frac{4}{27} \end{bmatrix}\]
Sendo a estatística \(T^2\) de Hotelling dada por
\[\begin{eqnarray*}T^2 &=& n(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_{0}}) \\ &=& 3 \times \begin{bmatrix} 8-9 & 6-5 \end{bmatrix}\begin{bmatrix} \frac{1}{3} & \frac{1}{9} \\ \frac{1}{9} & \frac{4}{27} \end{bmatrix}\begin{bmatrix} 8-9 \\ 6-5 \end{bmatrix}\\ &=& 3 \times \begin{bmatrix} -1 & 1 \end{bmatrix}\begin{bmatrix} -\frac{2}{9} \\ \frac{1}{27} \end{bmatrix} = \dfrac{7}{9}\end{eqnarray*}\]
Antes da amostra ser selecionada, a distribuição de \(T^2\) é dada por
\[\dfrac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha) = \frac{(3-1) \times 2}{(3-2)}F_{2,3-2} = 4F_{2,1}\]
Sendo assim, para testarmos as hipóteses
\[H_{0}:\mathbf{\mu} = \begin{bmatrix} 9 \\ 5 \end{bmatrix} \text{ vs. } H_{a} :\mathbf{\mu} \neq \begin{bmatrix} 9 \\ 5 \end{bmatrix}\]
Basta comparar o valor observado de \(T^2\) com o valor crítico da distribuição \(F_{2,1}\), ao nível de significância de \(\alpha = 0,05\), dado por \(4F_{2,1}(0,05) = 4 \times 199,5 = 798\)
\[T^2 = \dfrac{7}{9} = 0,778 << 798 = 4F_{2,1}(0,05)\] Logo, não existem evidências significativas a favor da rejeição de \(H_0\), ao nível de 5% de significância.
Exemplo 03: A transpiração de 20 mulheres saudáveis foi analisada. As observações são trivariadas, a saber: \(X_1\) - taxa de suor, \(X_2\) - conteúdo de sódio e \(X_3\) - conteúdo de potássio. Deseja-se testar, ao nível de significância de 10%, as seguintes hipóteses:
\[H_0: \mathbf{\mu} = \left[ \begin{array}{c} 4 \\ 50 \\ 10\\ \end{array} \right] \text{ vs } H_a: \mathbf{\mu} \neq \left[ \begin{array}{c} 4 \\ 50 \\ 10\\ \end{array} \right]\]
dados = read.table("https://raw.githubusercontent.com/tiagomartin/est022/refs/heads/main/dados/suor.dat", col.names = c("suor","sodio","potassio"))
dados %>% str()'data.frame': 20 obs. of 3 variables:
$ suor : num 3.7 5.7 3.8 3.2 3.1 4.6 2.4 7.2 6.7 5.4 ...
$ sodio : num 48.5 65.1 47.2 53.2 55.5 36.1 24.8 33.1 47.4 54.1 ...
$ potassio: num 9.3 8 10.9 12 9.7 7.9 14 7.6 8.5 11.3 ...
# Teste de normalidade multivariada
# H_0: dados seguem distribuição normal trivariada
mvshapiro_test(as.matrix(dados))
Generalized Shapiro-Wilk test for Multivariate Normality
data: as.matrix(dados)
MVW = 0.94446, p-value = 0.2567
p-valor = 0,2567 > 0,10 = \(\alpha\) → \(NRH_0\)
One Sample Hotelling T Square Test
Hotelling T Sqaure Statistic = 9.738773
F value = 2.905 , df1 = 3 , df2 = 17 , p-value: 0.0649
Descriptive Statistics
suor sodio potassio
N 20.00000 20.00000 20.000000
Means 4.64000 45.40000 9.965000
Sd 1.69687 14.13465 1.904641
Detection important variable(s)
Lower Upper Mu0 Important Variables?
suor 3.555292 5.724708 4 FALSE
sodio 36.364555 54.435445 50 FALSE
potassio 8.747476 11.182524 10 FALSE
p-valor = 0,0649 < 0,10 = \(\alpha\) → \(RH_0\)
\[\dfrac{n-p}{(n-1)p}T^2 = F_{p, n-p}\]
A rejeição da hipótese nula nos conduz ao problema de identificar quais variáveis são responsáveis pela rejeição dessa hipótese.
A região de confiança (RC) é uma generalização dos intervalos de confiança univariados.
A definição dessa região nos auxilia a identificar os componentes do vetor de parâmetros que são responsáveis pela rejeição de \(H_0\).
Seja \(\mathbf{\theta}\) um vetor de parâmetros populacional desconhecidos e \(\mathbf{\Theta}\) o conjunto de todos os possíveis valores de \(\mathbf{\theta}\). A RC é uma região de possíveis valores de \(\mathbf{\theta}\).
A RC é chamada de região com \(100(1 - \alpha)\%\) de confiança se, antes da amostra ser selecionada,
\[P(\mathbf{\theta} \in RC) = 1 - \alpha\]
Esta probabilidade é calculada sob o verdadeiro, porém desconhecido, valor de \(\mathbf{\theta}\).
A RC para o vetor de médias populacional \(\mathbf{\mu}\), quando se dispõe de uma amostra aleatória da distribuição \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\) é dada por,
\[P\left[ n(\bar{\mathbf{x}}-{\mathbf{\mu}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{X}}-{\mathbf{\mu}}) \leqslant \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha) \right] = 1 - \alpha\]
de forma que,
\[RC = \left\lbrace \mathbf{\mu} \in \mathbf{R}^p | n(\bar{\mathbf{x}}-{\mathbf{\mu}})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}}) \leqslant \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha) \right\rbrace \]
em que \(F_{p,n-p} (\alpha)\) denota o \(100(1 − \alpha)\)-ésimo percentil superior de uma distribuição \(F[p; (n−p)]\).
Temos que a região de confiança é dada pelo hiper elipsoide de eixos determinados pelos autovetores da matriz de covariância amostral \(\mathbf{S}\) e cujas medidas são proporcionais às raízes quadradas dos respectivos autovalores.
Para verificar se um dado vetor \(\mathbf{\mu}_0\) pertence à região de confiança, basta calcular o valor de
\[n(\bar{\mathbf{x}}-{\mathbf{\mu}_0})^{t}\mathbf{S}^{-1}(\bar{\mathbf{x}}-{\mathbf{\mu}_0})\]
e compará-lo com o valor de
\[\frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha)\]
Para \(p \geqslant 4\) não é possível representar visualmente a RC. No entanto, podemos calcular as medidas dos eixos do hiper elipsoide de confiança centrado em \(\mathbf{\bar x}\):
\[n(\mathbf{\bar x}-{\mathbf{\mu}})^{t}\mathbf{S}^{-1}(\mathbf{\bar x}-{\mathbf{\mu}}) \leqslant c^2 = \frac{(n-1)p}{(n-p)}F_{p,n-p}(\alpha)\]
Os semi-eixos têm medida
\[\sqrt{\lambda_i} \displaystyle{\frac{c}{\sqrt{n}}} = \sqrt{\lambda_i} \sqrt{\displaystyle{\frac{(n-1)p}{n(n-p)}F_{p,n-p}(\alpha)}}\]
unidades na direção de \(\mathbf{e}_i\), em que \(\lambda_i\) e \(\mathbf{e}_i\) são os autovalores e autovetores de \(\mathbf{S}\).
Começando no centro \(\mathbf{\bar x}\), os eixos do elipsoide de confiança são
\[\pm \left( \sqrt{\lambda_i} \sqrt{\displaystyle{\frac{(n-1)p}{n(n-p)}F_{p,n-p}(\alpha)}}\right) \mathbf{e}_i\]
em que \(\lambda_i\) e \(\mathbf{e}_i\) são os autovalores e autovetores de \(\mathbf{S}\), \(i = 1, \cdots, p\).
Exemplo 04: Uma empresa quer verificar se a média do desempenho dos funcionários em dois critérios-chave está de acordo com o padrão estabelecido pela diretoria.
Os critérios avaliados são:
A diretoria espera que a média populacional dos funcionários seja:
\[\mathbf{\mu}_0 = \left[ \begin{array}{c} 80 \\ 7,5 \end{array} \right]\]
dados = read.csv("https://raw.githubusercontent.com/tiagomartin/est022/refs/heads/main/dados/dados_desempenho.csv", row.names = 1)
dados %>% str()'data.frame': 20 obs. of 2 variables:
$ Produtividade: num 90.5 87.3 69.3 84.2 83.7 ...
$ Satisfacao : num 7.9 6.88 6.15 7.09 6.95 6.62 5.38 7.38 8.07 5.75 ...
One Sample Hotelling T Square Test
Hotelling T Sqaure Statistic = 67.55705
F value = 32.001 , df1 = 2 , df2 = 18 , p-value: 1.18e-06
Descriptive Statistics
Produtividade Satisfacao
N 20.000000 20.000000
Means 83.580500 6.431500
Sd 9.732401 1.026195
Detection important variable(s)
Lower Upper Mu0 Important Variables?
Produtividade 77.619031 89.541969 80.0 FALSE
Satisfacao 5.802916 7.060084 7.5 *TRUE*
p-valor = 0,00000118 < 0,05 = \(\alpha\) → \(RH_0\)
🤔 Qual variável foi a responsável?
Se a média esperada \(\mathbf{\mu}_0\) estiver fora da elipse, há indícios de que o desempenho real dos funcionários não está de acordo com a expectativa.
x_barra = colMeans(dados)
n = nrow(dados)
p = ncol(dados)
S = var(dados)/n
alpha = 0.05
c2 = (n-1)*p*qf(1-alpha,p,n-p)/(n-p)
ci.prod=ci_mean(dados$Produtividade, type = "t")$interval
ci.satis=ci_mean(dados$Satisfacao, type = "t")$interval
# Cálculo da elipse de confiança (nível 95%)
Radius = sqrt(eigen(S)$values) * sqrt(c2)
segments = 51
angles = (0:segments) * 2 * pi/segments
unit.circle = cbind(cos(angles), sin(angles))
Q = sweep(eigen(S)$vectors, 2, Radius, "*")
elipse_conf = data.frame(sweep(unit.circle %*% t(Q), 2, as.vector(x_barra), "+"))
names(elipse_conf) = c("Produtividade", "Satisfacao")
# Criando o gráfico
dados %>%
ggplot() +
geom_point(aes(x = Produtividade, y = Satisfacao), color = "blue", alpha = 0.6) +
geom_path(data = elipse_conf, aes(x = Produtividade, y = Satisfacao), color = "red") +
geom_point(aes(x = x_barra[1], y = x_barra[2]), color = "red", size = 3) +
geom_point(aes(x = mu0[1], y = mu0[2]), color = "black", size = 3, shape = 17) +
geom_hline(yintercept=ci.satis, linetype="dashed", color = "red") +
geom_vline(xintercept = ci.prod,linetype="dashed", color = "red")+
labs(title = "Elipse de Confiança 95% para a Média Populacional",
x = "Produtividade (tarefas/dia)",
y = "Satisfação no Trabalho (escala 1-10)") +
theme_minimal()Quando \(p \geqslant 4\), não conseguimos visualizar graficamente a região de confiança, de forma que necessitaremos de um procedimento alternativo para descobrir qual componente do vetor de parâmetros foi o responsável pela rejeição de \(H_0\).
A ideia básica é fazermos intervalos de confiança para cada componente, ou para combinações lineares dos componentes do vetor de parâmetros.
Sejam \(\mathbf{x} \sim N_p(\mathbf{\mu}, \mathbf{\Sigma})\) e a combinação linear
\[Z = \mathbf{l}^t \mathbf{x} = l_1X_1 + l_2X_2 + \cdots + l_pX_p\]
Temos que
\[\mu_Z = E(Z) = \mathbf{l}^t \mathbf{\mu} \hspace{0.5cm} \text{e} \hspace{0.5cm} \sigma^2_Z = \mathbf{l}^t \mathbf{\Sigma} \mathbf{l}\]
de forma que
\[Z \sim N(\mathbf{l}^t \mathbf{\mu}, \mathbf{l}^t \mathbf{\Sigma} \mathbf{l})\]
Logo, se \(\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\) é uma amostra aleatória de uma população \(N_p(\mathbf{\mu}, \mathbf{\Sigma})\), segue que \(Z_1, Z_2, \cdots, Z_n\), definidos por
\[Z_i = \mathbf{l}^t \mathbf{x}_i, \hspace{0.5cm} i = 1, 2, \cdots n\]
é uma amostra aleatória de uma população \(N(\mathbf{l}^t \mathbf{\mu}, \mathbf{l}^t \mathbf{\Sigma} \mathbf{l})\).
Da teoria normal univariada, para \(\mathbf{l}\) fixo e \(\sigma^2\) desconhecido, temos que um intervalo de \(100(1-\alpha)\%\) de confiança para \(\mu_Z = \mathbf{l}^t \mathbf{\mu}\) é baseado na razão t-Student, isto é,
\[ t = \displaystyle{\frac{\bar{z} - \mu_Z}{s_Z/\sqrt{n}}} = \displaystyle{\frac{\sqrt{n}(\mathbf{l}^t \bar{\mathbf{x}} - \mathbf{l}^t \mathbf{\mu})}{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}} \]
O que conduz a,
\[\bar{z} - t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{s_Z}}{\sqrt{n}}} \leqslant \mu_Z \leqslant \bar{z} + t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{s_Z}}{\sqrt{n}}}\]
ou,
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
Podemos construir vários intervalos de confiança sobre combinações lineares dos componentes do vetor \(\mathbf{\mu}\), cada um associado com um coeficiente de confiança \(1 − \alpha\), escolhendo diferentes vetores de constantes \(\mathbf{l}\). Porém, o coeficiente de confiança conjunto do grupo de intervalos resultantes não seria mais \(1 − \alpha\).
É desejável então, associar um coeficiente de confiança coletivo de \(1 − \alpha\) aos intervalos de confiança que podem ser gerados para diferentes escolhas de \(\mathbf{l}\) (Bônus).
Ônus: intervalos que têm amplitudes maiores (mais largos, menos precisos) do que os intervalos apresentados anteriormente via a distribuição amostral t-Student com \(n − 1\) graus de liberdade.
Uma solução para o problema de cobertura é a construção dos intervalos simultâneos. Dado um conjunto de dados observados \(\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\) e uma particular escolha do vetor \(\mathbf{l}\), temos que o intervalo de confiança
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + t_{\left[ (n-1);\frac{\alpha}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
representa o conjunto de valores de \(\mathbf{l}^t \mathbf{\mu}\) para os quais
\[|t| = \left|\displaystyle{\frac{\sqrt{n}(\mathbf{l}^t \bar{\mathbf{x}} - \mathbf{l}^t \mathbf{\mu})}{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}}\right| \leqslant t_{\left[ (n-1);\frac{\alpha}{2}\right]}\]
Ou, equivalentemente,
\[t^2 = \displaystyle{\frac{n(\mathbf{l}^t \bar{\mathbf{x}} - \mathbf{l}^t \mathbf{\mu})^2}{\mathbf{l}^t \mathbf{S} \mathbf{l}}} = \displaystyle{\frac{n\left[ \mathbf{l}^t( \bar{\mathbf{x}} - \mathbf{\mu})\right] ^2}{\mathbf{l}^t \mathbf{S} \mathbf{l}}} \leqslant t^2_{\left[ (n-1);\frac{\alpha}{2}\right]}\]
Uma região de confiança simultânea é dada pelos valores de \(\mathbf{l}^t \mathbf{\mu}\) tais que \(t^2\) é relativamente pequeno para todas as escolhas de \(\mathbf{l}\).
Parece razoável esperar que a constante \(t_{\left[ (n-1);\frac{\alpha}{2}\right]}\) seja trocada por um valor maior, \(c^2\), quando as declarações forem desenvolvidas para muitas escolhas de \(\mathbf{l}\).
A constante \(c^{2}\) deve ser obtida para manter o valor global do coeficiente de confiança \(1-\alpha\) simultâneo para todas as escolhas de \(\mathbf{l}\).
Pode ser demonstrado que, tomando
\[c^{2}= \displaystyle{\frac{(n-1)p}{n-p}F_{p,n-p}(\alpha)}\]
tem-se a garantia do valor global para o coeficiente de confiança igual a \(100(1-\alpha)\%\) para qualquer escolha não-nula do vetor \(\mathbf{l}\).
É conveniente nos referirmos aos intervalos simultâneos como intervalos \(T^2\), uma vez que a probabilidade de cobertura é determinada pela distribuição de \(T^2\)
Assim, os intervalos de confiança simultâneos (ou intervalos \(T^2\)) são dados por:
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - \sqrt{\displaystyle{\frac{(n-1)p}{n-p}F_{p,n-p}(\alpha)}} \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + \sqrt{\displaystyle{\frac{(n-1)p}{n-p}F_{p,n-p}(\alpha)}} \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
Exemplo 05: Voltando ao exemplo anterior referente à uma empresa que quer verificar se a média do desempenho dos funcionários em dois critérios-chave (produtividade e satisfação) está de acordo com o padrão estabelecido pela diretoria.
Encontre os intervalos de confiança simultâneos \(T^2\) de Hotelling e compare com os intervalos univariados.
lim = cbind(x_barra - sqrt(c2)*sqrt(diag(S)),
x_barra + sqrt(c2)*sqrt(diag(S)))
# Criando o gráfico
dados %>%
ggplot() +
geom_point(aes(x = Produtividade, y = Satisfacao), color = "blue", alpha = 0.6) +
geom_path(data = elipse_conf, aes(x = Produtividade, y = Satisfacao), color = "red") +
geom_point(aes(x = x_barra[1], y = x_barra[2]), color = "red", size = 3) +
geom_point(aes(x = mu0[1], y = mu0[2]), color = "black", size = 3, shape = 17) +
geom_hline(yintercept = ci.satis, linetype="dashed", color = "red") +
geom_vline(xintercept = ci.prod,linetype="dashed", color = "red")+
geom_hline(yintercept = lim[2,], linetype="dashed", color = "blue") +
geom_vline(xintercept = lim[1,], linetype="dashed", color = "blue")+
labs(title = "Elipse de Confiança 95% para a Média Populacional",
x = "Produtividade (tarefas/dia)",
y = "Satisfação no Trabalho (escala 1-10)") +
theme_minimal()Uma alternativa aos intervalos simultâneos são os Intervalos de Confiança com Correção de Bonferroni, baseada na Desigualdade de Bonferroni.
Sejam \(A_1, A_2, \cdots, A_m\) uma coleção de eventos em um espaço de probabilidade tais que \(P(A_i) = 1 - \alpha_i\), \(i = 1, \cdots, m\). Então,
\[ \begin{eqnarray*} P(\cap_{i = 1}^m A_i) &=& 1 - P \left(``\text{pelo menos um } A_i \text{ não ocorre}"\right) \\&\geqslant& 1 - \displaystyle{\sum_{i = 1}^m P(A_i^c)} = 1 - \displaystyle{\sum_{i = 1}^m \alpha_i} \\&=& 1 - (\alpha_1 + \alpha_2 + \cdots + \alpha_m) \end{eqnarray*} \]
A desigualdade de Bonferroni nos permite o controle da taxa de erro global \(\alpha_1 + \alpha_2 + \cdots + \alpha_m\), independente da estrutura de correlação que está por trás dos intervalos de confiança.
Assim, se o problema envolve a construção de \(m\) intervalos de confiança, a ideia é fazer
\[\alpha^{\ast}=\dfrac{\alpha}{m}\]
e tomar os intervalos separadamente, dados por
\[ \begin{eqnarray} \mathbf{l}^t \bar{\mathbf{x}} - t_{\left[ (n-1);\frac{\alpha^{\ast}}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \leqslant \mathbf{l}^t \mathbf{\mu} \leqslant \mathbf{l}^t \bar{\mathbf{x}} + t_{\left[ (n-1);\frac{\alpha^{\ast}}{2}\right] } \displaystyle{\frac{\sqrt{\mathbf{l}^t \mathbf{S} \mathbf{l}}}{\sqrt{n}}} \end{eqnarray} \]
Observe agora, que o coeficiente de confiança coletivo vale pelo menos,
\[1 - \left( \underbrace{\displaystyle{\frac{\alpha}{m}} + \displaystyle{\frac{\alpha}{m}} + \cdots + \displaystyle{\frac{\alpha}{m}}}_{\text{m termos}}\right) = 1 - \alpha\]
Estes intervalos podem, de forma mais apropriada, serem comparados com os intervalos \(T^2\).
Exemplo 06: Voltando ao exemplo anterior referente à uma empresa que quer verificar se a média do desempenho dos funcionários em dois critérios-chave (produtividade e satisfação) está de acordo com o padrão estabelecido pela diretoria.
Encontre os intervalos de confiança com correção de Bonferroni e os compare com os intervalos simultâneos \(T^2\) de Hotelling e com os intervalos univariados.
m=2
alpha_corr = alpha/m
ci_bonf.prod=ci_mean(dados$Produtividade, type = "t", probs=c(alpha_corr/2, 1-alpha_corr/2))$interval
ci_bonf.satis=ci_mean(dados$Satisfacao, type = "t", probs=c(alpha_corr/2, 1-alpha_corr/2))$interval
# Criando o gráfico
dados %>%
ggplot() +
geom_point(aes(x = Produtividade, y = Satisfacao), color = "blue", alpha = 0.6) +
geom_path(data = elipse_conf, aes(x = Produtividade, y = Satisfacao), color = "red") +
geom_point(aes(x = x_barra[1], y = x_barra[2]), color = "red", size = 3) +
geom_point(aes(x = mu0[1], y = mu0[2]), color = "black", size = 3, shape = 17) +
geom_hline(yintercept = ci.satis, linetype="dashed", color = "red") +
geom_vline(xintercept = ci.prod,linetype="dashed", color = "red")+
geom_hline(yintercept = lim[2,], linetype="dashed", color = "blue") +
geom_vline(xintercept = lim[1,], linetype="dashed", color = "blue") +
geom_hline(yintercept = ci_bonf.satis, linetype="dashed", color = "#006400") +
geom_vline(xintercept = ci_bonf.prod,linetype="dashed", color = "#006400")+
labs(title = "Elipse de Confiança 95% para a Média Populacional",
x = "Produtividade (tarefas/dia)",
y = "Satisfação no Trabalho (escala 1-10)") +
theme_minimal()Novamente, vamos revisar o caso univariado. Sejam uma amostra aleatória \(x_{11}, x_{12}, \cdots, x_{1n_1}\) de uma população \(N(\mu_1, \sigma^2_1)\) e uma segunda amostra \(x_{21}, x_{22}, \cdots, x_{2n_2}\) de uma população \(N(\mu_2, \sigma^2_2)\). Vamos assumir que as amostras são independentes e que \(\sigma^2_1 = \sigma_2^2 = \sigma^2\), com \(\sigma^2\) desconhecido.
Sejam as seguintes hipóteses:
\[H_0: \mu_1 = \mu_2 \text{ vs. } H_a: \mu_1 \neq \mu_2\]
Uma estatística de teste apropriada para este tipo de situação é
\[t = \dfrac{\bar{x_1}- \bar{x}_2}{s_p\sqrt{\dfrac{1}{{n_1}}+\dfrac{1}{n_2}}} \stackrel{H_0}{\sim} t_{n_1 + n_2 - 2}\]
em que
\[s_p^2 = \dfrac{(n_1-1)s_1^2 + (n_2-1)s_2^2}{n_1 + n_2 - 2}\]
Considere agora uma amostra aleatória \(\mathbf{x_1} = \{X_{11}, X_{12}, \cdots, X_{1n_1}\}\) de uma população \(N_p(\mathbf{\mu_1}, \mathbf{\Sigma_1})\) e uma segunda amostra \(\mathbf{x_2} = \{X_{21}, X_{22}, \cdots, X_{2n_2}\}\) de uma população \(N_p(\mathbf{\mu_2}, \mathbf{\Sigma_2})\)
Assumiremos que as duas amostras são independentes e que \(\mathbf{\Sigma_1} = \mathbf{\Sigma_2} = \mathbf{\Sigma}\).
Queremos testar as seguintes hipóteses:
\[H_0: \mathbf{\mu}_1 = \mathbf{\mu}_2 \text{ vs. } H_a: \mathbf{\mu}_1 \neq \mathbf{\mu}_2\]
Quando as matrizes de covariâncias populacionais são iguais, as hipóteses podem ser testadas utilizando a seguinte estatística,
\[T^2 = (\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)^t \left[ \left( \dfrac{1}{n_1} + \dfrac{1}{n_2} \right) \mathbf{S}_p \right]^{-1}(\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)\]
em que
\[\mathbf{S}_p = \dfrac{(n_1-1)\mathbf{S}_1 + (n_2-1)\mathbf{S}_2}{n_1 + n_2 - 2}\]
Sob \(H_0\), temos que
\[T^2 \stackrel{H_0}{\sim} \displaystyle{\frac{p(n_1 + n_2 - 2)}{n_1 + n_2 - p - 1}} F_{(p,n_1+n_2-p-1)}\]
Ao nível \(\alpha\) de significância, rejeita-se \(H_{0}\) em favor de \(H_{a}\) se
\[T^{2}\geq \displaystyle{\frac{p(n_1 + n_2 - 2)}{n_1 + n_2 - p - 1}} F_{(p,n_1+n_2-p-1)}(\alpha)\]
em que \(F_{(p,n_1+n_2-p-1)}(\alpha)\) denota o \(100(1 − \alpha)\)-ésimo percentil superior de uma distribuição \(F[p,(n_1+n_2-p-1)]\).
\[\dfrac{n_1+n_2-p-1}{(n_1+n_2-2)p}T^2 = F_{p, n_1+n_2-p-1}\]
Exemplo 06: Quatro testes psicológicos foram aplicados em 32 homens e 32 mulheres. As variáveis mensuradas foram:
Teste as hipóteses de que
\[H_0: \mathbf{\mu}_{\text{Homens}} = \mathbf{\mu}_{\text{Mulheres}} \text{ vs. } H_a: \mathbf{\mu}_{\text{Homens}} \neq \mathbf{\mu}_{\text{Mulheres}}\]
dados = read.csv("https://raw.githubusercontent.com/tiagomartin/est022/refs/heads/main/dados/testes_psico.csv", sep = ";") %>% mutate(Sexo=ifelse(Sexo=="M", 1,2))
dados %>% str()'data.frame': 64 obs. of 5 variables:
$ X1 : int 15 17 15 13 20 15 15 13 14 17 ...
$ X2 : int 17 15 14 12 17 21 13 5 7 15 ...
$ X3 : int 24 32 29 10 26 26 26 22 30 30 ...
$ X4 : int 14 26 23 16 28 21 22 22 17 27 ...
$ Sexo: num 1 1 1 1 1 1 1 1 1 1 ...
res = TwoSamplesHT2(dados %>% dplyr::select(-Sexo), group=dados$Sexo, alpha = 0.05, Homogenity = TRUE)
summary(res) Two Independent Samples Hotelling T Square Test
Hotelling T Sqaure Statistic = 97.6015
F value = 23.22 , df1 = 4 , df2 = 59 , p-value: 1.46e-11
Descriptive Statistics (The First Group)
X1 X2 X3 X4
Means 15.968750 15.906250 27.187500 22.750000
Sd 2.278715 3.631043 5.354754 4.079848
Descriptive Statistics (The Second Group)
X1 X2 X3 X4
Means 12.343750 13.906250 16.656250 21.937500
Sd 3.022596 4.313216 5.480813 5.291122
Detection important variable(s)
Lower Upper Important Variables?
X1 1.443739 5.806261 *TRUE*
X2 -1.248920 5.248920 FALSE
X3 6.115836 14.946664 *TRUE*
X4 -3.037608 4.662608 FALSE
p-valor = 1,46e-11 < 0,05 = \(\alpha\) → \(RH_0\)
Para o caso em que \(\mathbf{\Sigma}_1 \neq \mathbf{\Sigma}_2\), ou seja, as matrizes de covariâncias populacionais se diferem, sugere-se o uso da estatística
\[T^2 = (\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)^t \left[ \left( \dfrac{\mathbf{S}_1}{n_1} + \dfrac{\mathbf{S}_2}{n_2} \right) \right]^{-1}(\mathbf{\bar{x}}_1 - \mathbf{\bar{x}}_2)\]
Para grandes amostras, pode-se utilizar
\[T^2 \stackrel{H_0}{\sim} \chi^2_p,\] sendo \(\chi^2_p\) a distribuição qui-quadrado com \(p\) graus de liberdade.
A hipótese nula deve ser rejeitada se \(T^2 \geq \chi ^2_p(\alpha)\) para o nível de significância \(\alpha\).